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Abstract 



Shape estimation and object reconstruction are common problems in image analysis. Mathematically, 
' viewing objects in the image plane as random sets reduces the problem of shape estimation to inference 

C/3 . about sets. Currently existing definitions of the expected set rely on different criteria to construct the 

expectation. This paper introduces new definitions of the expected set and the expected boundary, 
^ based on oriented distance functions. The proposed expectations have a number of attractive properties, 

\^ ' including inclusion relations, convexity preservation and equivariance with respect to rigid motions. The 

\l ■ paper introduces a special class of separable oriented distance functions for parametric sets and gives 

the definition and properties of separable random closed sets. Further, the definitions of the empirical 
^ mean set and the empirical mean boundary are proposed and empirical evidence of the consistency of the 

' boundary estimator is presented. In addition, the paper gives loss functions for set inference in frequentist 

framework and shows how some of the existing expectations arise naturally as optimal estimators. The 
proposed definitions of the set and boundary expectations are illustrated on theoretical examples and 
real data. 



1 Introduction 

Boundary reconstruction and shape estimation are frequently encountered problems in image analysis. For 
example, it is often of interest to determine a characteristic shape of a cell, reconstruct tissue boundary in 
medical images, or estimate a probable area in the earthquake disaster zone. Mathematically, the objects 
of interest can be viewed as sets, whilst inherent stochasticity of the acquisition process turns them into 
random entities. The problems of boundary reconstruction and shape estimation thus reduce to inference 
about random sets. 



Early results in the the ory of random se ts date back to IChoauetl (|l953f ) with a thorough mathemat 



( Molchanov 



leal treatment app earing in iMatheronI (|l975l ) and the most recent developments in the field presented in 



2005( 1. Because the family of closed sets is nonlinear, there is no natural way to define the 
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expected set. Currently, there e xist a number of de finitions of the expectation, though neither can be spoken 
of as the best. As pointed out in lMolchanovl ([20051), the definition of the expectation depends on set features 
that are important to emphasize. 

Many of the existing definitions are not based on a random set directly, but use an embedding into a 
function space. For example, the Vorob'ev definition uses the representation of a set given by its characteristi c 
function to construct an expectation that is optimal with respect to Lebesgue measure (IVorob'evl . 119841 ) . 
The distance-average approach maps a set into the space of distance functions, so that the expected set is 
optimal among all of the l evel s ets of the expected distance function with respect to a predefined metric 
(jBaddelev and Molchanovi llQQsl) . 

Currently existing definitions of the expected set do not address the question of boundary estimation, 
which is often of primary interest in practice. In this paper, we introduce a new definition of the expected 
set and the expected boundary based on oriented distance functions (ODFs). As opposed to a conventional 
distance function, an O DF takes into account the set and its complement, providing a more informative 
representation of a set (jPelfour and Zolesiol 120011 ). In comparison with currently existing definitions, the 
new expectations have a number of attractive theoretical properties. In particular, they satisfy inclusion 
relations, preserve convexity and remain equivariant with respect to rigid motions. We introduce a new 
class of ODFs for parametric sets and show the connection between the expectation of a parametric random 
closed set and the expected parameter. We also outline a general framework for set inference and show how 
different expectations arise as natural estimators for different loss functions. 



2 Random Closed Sets 

2.1 Notation 

For set A C K*^, denote by A'^, dA, cl A, int A, X{A) its complement, boundary, closure, interior, and Lebesgue 
measure, respectively. We write x = {xi, . . . ,Xd) for point x G R'' and use |- | for the standard Euclidean 
norm. Let T (resp. K.) be the family of closed (resp. compact) subsets of R"* and let the tripple (il, £/, P) 
denote the probability space. 

2.2 Definitions and Examples 

Definition 2.1. A random closed set is the mapping A : il i-^ T such that for every compact set K G JC, 

{uj : A{uj) n K 7^ 0} G ^. 

Here, we consider random closed sets taking values in R''. Throughout this paper, we write "r.c.s." for 
a random closed set, though some texts prefer RACS. We write A for an r.c.s. A(a;) and reserve capital 
roman letters for closed subsets of R''. 

Random closed sets give rise to several random variables and random functions, including the character- 
istic function 

, , f 1, if a; G ^, 
[0, lix ^ A, 
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the Lebesgue measure X{A) = J^^ X^(x)dx, and the distance function from point x G M.'^ to A, 



dA{x) 



+00, A = %, 



(2) 



where p{x,y) is a metric on M''. Some examples of p include the L2 Euclidean metric, the Lqo Chebyshev 
metric, and the Li Manhattan metric. Throughout this paper, we use the Euclidean distance function with 

pix,y) = \x-y\ in 

2.3 Expectations of Random Closed Sets 

Because the space T of all closed subsets of R'' is nonlinear, there is no natural way to define the expected 
set. In this section, we briefly review some of the existing definitions and refer for more details to ( Aumannl . 



1965 



Artstein and Vitak 



1975 



Vorob'ev 



1984 



Baddelev and Molchanov 



1998 



Molchanov 



20051) 



Definition 2.2 (Selection expectation). A random element ^ is called the selection of A, if it belongs to 
A with probability one. A selection ^ is integrable, if E |^| < 00. The selection expectation. Eg [A], is the 
closure of the set of the expectations of all integrable selections of A, 

Es[A] = cl {E[^] : ^ is a selection and E |^| < 00}. 



The selection expectation depends on the structure of the probability space (jMolchanovl . |2005| . see Example 
1.14, Section 2). In addition, if the probability space is nonatomic. Eg [A] is necessarily convex, even for 
nonconvex deterministic sets. 

In the conte xt of image analysis, the Vorob'ev expectation is perhaps the most intuitive construction 



( Vorob'ev 



19841 ). Consider an r.c.s. A C P. For every point x (zV, we assign a probability mass depending 



on whether or not x G A. For example, assigning one to all points in A and zero to points in A"^ converts the 
observed set into a binary image. Averaging binary images over all realizations of A, we obtain a gray-scale 
image, where the intensity at each pixel is the probability of the point being in A. The gray-scale image is 
not a binary image, unless A is deterministic. The Vorob'ev definition then provides a criterion to construct 
an expected set that is optimal with respect to Lebesgue measure. 

More precisely, given the characteristic function Xa of an r.c.s. A, define the coverage function by 



E[XA{x)]^P{xeA). 
The excursions sets of the coverage function arc given by 

Au = {xeM.'^ : E[XAx)] >u}, q€ [0, 1]. 



(3) 



(4) 



Definition 2.3 (Vorob'ev expectation). The expectation, Ev[A], of a random closed set A is the excursion 
set Ag, where q G [0, 1] is such that 



A(A„) < E[A(A)] < A(AJ for all w > q. 



(5) 
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Note that the optimality criterion ([5]) ensures that the Vorob'ev expectation ignores sets of measure zero. 

The distance-average expecta tion is based on a representation of a set given by some function /a : 25 i-^ M 
(jBaddelev and Molchanovl 119981 ) . We call /a a representative function of the set A. Examples of / include 
the distance function, the ODF, the characteristic function and others. For an r.c.s. A with representative 
function /a, let E[/a(x)] be the expected value of /a at x, assuming it exists. 

Definition 2.4 (Distance-average expectation). For a compact window W, define a (pseudo-) metric 
m(A, B) ~ tnw(/A(-)' /b(-))- The distance-average expectation, Eoa[A], of an r.c.s. A is the level set 



A„ = {x e W: E[/A(a;)] < u}, where u = arginf m(E[/A], /aJ 



(6) 



Examples of pseudo metric m in (|6l) include Lg distances and their variates, Baddeley's A^-distancc, or 
any custom definitions IjBaddelevi . 119921 ). Note that the pseudometric m is computed over the window W, 
which can be the entire domain V or its subset. For brevity, we omit the subscript W whenever possible. 

The distance-average expectation strongly depends on the choice of the repre sentative function /, window 
W, (pseudo-) metric m and parameters of m (jBaddelev and Molchanovl Il998l ). The Vorob'ev expectation 
ca n be viewed as a special case of the distance-average definition with fA.{x) = 1 — xa{x) and distance 
m (jBaddelev and Molchanovl . 1 19981 see Example 5.14). 

The linearization approach gives a unified view of the existing expectations (jMolchanovl 120051 ). In 
particular, consider a mapping f : J- t-^ y, where y is some Banach space. Let /a be the image of A in 
y and assume that E[|/a(.t)|] < oo. If there exists a unique F G T such that /i?(x) = E[/A(a;)], then 
declare E[A] = F. In general, the inverse image rarely exists, and we need to define the criterion so that the 
resulting expectation is optimal in some sense. This is done as follows. Let J) be a pseudometric in y, then 
the expectation of an r.c.s. A is given by 



E[A] = argmin5(E[C(A)],C(F)). 



(7) 



Minimizing over the entire family is often nontrivial. Therefore, we can consider optimizing Z)(-, •) over 
a subfamily H C T oi candidate sets. Clearly, the resulting expectation depends on the choice of Banach 
space 3^, mapping /, metric 0, subfamily Ti and any implicit parameters that enter into the calculation. 

For bounded closed covex sets, the selection expectation is an example of a linearization approach, where 
the embedding into a Banach space is realised by support functions of random sets on the unit sphere 
and the expected set is given by the support function of the selection expectation. Another example is 
the Vorob'ev definition, which is based on the mapping given by Xa ^-nd gives the expected set Ev[A] 
which is optimal with respect to Lebesgue measure; the family of candidate sets Ti. consists of excursion 
sets (jlj of the coverage function. The distance-average expectation is based on the embedding into the 
space of representative function and the res ulting average is optimal among the level sets of t he expected 
representation with respect to metric m(-, •) (jBaddelev and Molchanovl . Il998t lMolchanovll2005l ). 



2.4 Oriented Distance Functions 

In comparison to a distance function ([2), an ODF takes into account both the set and its compliment, thus 
reflecting the interior and exterior properties. In addition, the ODF provides a level-set description of the 
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boundary. 

Definition 2.5. The ODF from point x to set A C M'^ with dA^%\s defined by 

hA{x) = dA{x) - dA-ix), for aU x e R''. (8) 
Using basic properties of the distance function, the expression ([5]) can be written as 

{dA{x) = ddA{x), X G int A'^ 
0, xedA 
—dA^{x) ~ -~doA{x), X e int A . 

Below we summarize some properties of the ODFs and refer for proofs and details to 
poo 11) . 

Proposition 2.1. Let A and B he some subsets of cl T) with dA, dB ^ 0. 

1. The ODF provides a level-set description of a set, i.e 

clA = {x:bA(x)<{)}, dA= {x:hA{x) ^Q}. 

2. A D B implies Ba < bs- 

3. bA <bB mV iff clB d clA and cl A" C cl B^ . 

4. bA = bs in V iff cl A ~ cl B and dA = dB. 

5. For a convex set A with dA 7^ 0, Ba = b^iA is a convex function in . 

6. The function bA{x) is uniformly Lipschitz in with constant one, i.e. 

\bA(.y) - bA{x)\ <\y - x\ for any x,y e R"^. 

In addition, bA{x) is (Frechet) differentiahle with \^bA{x)\ < 1 a.e. inW^. When it exists, the gradient 
of the ODF coincides with the outward unit normal to the boundary. 

Note that the ODF identifies the set up to its closure and boundary. In other words, ODF describes 
equivalence classes of sets with identical closure and boundary. 

Example 2.1. For a singleton A = {9}, 9 G M'', the ODF coincides with the distance function and 
= dA{x) = |a; — ^1 for all x e R'^. Note that bA^x) > for all x ^ 6 and the isocontours of bA{x) are 
spheres centered at 9; Figure [1] 1st plot. 

Example 2.2. Consider a closed ball in R'^ centered at the origin with radius 9, i.e. A ^ {x gM.'^ : \x\ < 9}. 
The distance function of A is ^^(2;) = |a;| — 9, for x G A, and zero, elsewhere. Hence, the ODF has a form 
bA{x) = I a; I —6 for all x G R''. The isocontours of bA{x) are spheres centered at the origin; Figure [H 2nd 
plot. 
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(9) 



Delfour and Zolesio 



Example 2.3. The distance function of a half plane set A = {x: xi < di} is given by dA{x) = xi — 9i for 
X G A"^. The ODF has a form xi ~ 9i for all x G with isocontours given by vertical lines parallel to 
the boundary; Figure [1] 3rd plot. 

Example 2.4. The ODF of an upper-plane set A ~ {x: X2 > xitan^} in is b^ix) — \x\ sin {9 — uj) , 
where to = arcsin (a;2|a;|^^); Figure [U 4th plot. 




Figure 1: From left to right: The ODF of a singleton {9}, a unit disc, a half-plane and an upper half-plane 
in M^; the zero-level isocontour (white) shows the boundary of the set. 



3 Expectations of Random Sets via ODFs 
3.1 Definition and Properties 

For an r.c.s. A with dA ^ 0, the ODF bAix) at point a; G R'' is a random variable taking values in R. 
Denote by E[6a(2;)] the expected value of bA^x), assuming it exists. We begin by proving some properties 
of the expected ODF. 

Proposition 3.1. Consider an r.c.s. A such that dA 7^ a.s. Then, 

1. The expected ODF is uniformly Lipshitz with constant 1. 

2. //E[|6a(2;o)|] < 00 for some xq G T>, then 'Ei\\bA{x)\\ < 00 for any x £T>. 

3- IfbA_{x) is differentiate at some xq G V a.s., then \7 Fi[bA^(xo)] exists and |V E[6A(a;o)]| < 1- 

Proof dH) Foranya:,y gR^ | E[feA (?/)]- E[6A(a;)]| < E[|6a(2/) - &a(x)|] <\y-x\. 

(21) Using the triangle inequality, boundedness of T> and Lipshiptz continuity of E[6A(a;)], we obtain 
E[|6A(a;)|] < \x - xo\ + E[|6A(a;o)|] < 00. 

([3]) Because VbA{xo) exists a.s. and |V6yi(xo)| < 1, we can interchange the differentiation on the space 
n and integration with respect to x, so the existence of VE[&A(a;o)] follows. From Vb_x{xo) < 1, we have 
|VE[&a(2:o)]| <E[|V6a(.to)|] < 1. □ 

Definition 3.1 (ODF expectation). Let A be a random closed set such that dA ^ a.s. and suppose that 
E[&A(a;o)] < 00 for some xq G T>. Then the expectation of an r.c.s. A is the zero-level set of the expected 
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ODF and expectation of the boundary of A is the zero-level isocontour of the expected ODF, i.e. 



E[A] = {x : E[6a(x)] < 0}, (10) 
E[dA] = {x:E[bA{x)]^Q}. (11) 

We refer to E[A] as "the expected set" or "the ODF expectation" of A. Note, that the expected set defined 
by (|3.ip can be an empty set. Later we give examples of such random sets and discuss the appropriateness 
of the outcome. 

Theorem 3.1 (Inclusion properties). Consider r.c.s. A,B C cl D with 9A, 9B 7^ and assume that 
E[bA_{x)] and E[6b(2:)] exist for some x clV. That is, suppose E[A] and E[B] are well-defined. Then, 

1. The expectation E[A] is a closed set. 

2. If A. ^ A a.s., for some deterministic A Cz then E[A] = A and E[9A] = dA. 

3. Suppose that an r.c.s. A satisfies dA. = B a.s. for some deterministic set B, then E[9A] D B. 

4. For an r.c.s. A, 9E[A] C E[c}A]. 

5- IfE[bA{x)] > E[bB{x)] a.s., then E[A] C E[B]. 

6. // A C B a.s., then E[A] C E[B]. 

Proof. (HI) The closedness of E[A] follows from the continuity of E[6a]- 
([2]) is straightforward. 

^ It suffices to note that for all points x g B, bsix) — b^ix) = 0, so that E[6A(a^)] = and B C E[9A]. 
(|4]) For any x e 9E[A], E[&a(-'i;)] = and hence x G E[9A]. Example 13.21 shows that the reverse relation 
does not hold. 

dH) For any x e E[A], > E[bA{x)] > E[bB{x)]. Thus x G E[B] and E[A] C E[B]. 

([S]) The proof follows immediately from the previous statement and Property [3] of the ODFs, which 
implies that for A C B, 6a (a;) > b-B{x) in V and hence, E[6a(-t)] > E[feB(a;)] a-s.. □ 

Corollary 3.1. If an r.c.s. A Z) U a.s. for some deterministic set U C cll), then E[A] D clU. Similarly, 
if an r.c.s. A CW a.s. for some deterministic set W C clV, then E[A] C cl W. 

Proof. Follows from statements ([2]) and ^ of Theorem l3.ll □ 



Although the inclusion properties of Theorem 13 . 1 1 and CoroU ar v 13 . 1 1 seein natu ral, they do not necessarily 
hold for other definitions, including the selection expectation ijMolchanovl . 120051 ) and the distance-average 
definition. The latter does not satisfy a weaker version of inclusion in the sense that the resulting expectation 
i s not necessarily contained in the convex hull of the union of all possible realizations of a random set 
(jBaddelev and Molchanovi 119981) . 

Example 3.1. (random singleton) Let ^ be a random variable and consider an r.c.s. A = {^}. Then, for 
any x eV, E[6A(a;)] = E[|a; - ^| > 0] and hence E[A] = 0. 
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The last result does not agree with suggestions in lMolchanovl lj2005r ). where E[{^}] = {E[^]} is listed as one 



of the desirable properties of the expected set. Indeed, the selection expectation s atisfies the suggested prop- 



erty a nd so does the distance-average definition, for certain choices of tn, / and W (jBaddelev and Molchanov 



19981 see Example 5.4). However, the ODF definition of the expected set was primarily motivated by prob- 
lems in image analysis, where random singletons are often regarded as noise. Example 13.11 highlights the 
natural denoising property of the ODF expectation in a sense that any speckles, unless observed a.s., are not 
present in the average reconstruction. We argue that if random singletons arc the sets of primary interests, 
they should be viewed as random vectors rather than random sets, in which case conventional definitions of 
the expectation and dispersion are applicable. 

More generally, Example l3.1l can be extended to random sets with zero Lebesgue measure. In particular, 
for a r.c.s. A with int A = and A(A) = 0, we have &A(a;) > for all x G A'^ and bA{x) = for all x S OA. 
Thus, unless x £ OA almost surely, E[6A(a;)] > and the expected set E[A] = {x: .t G A almost surely}. 

In the rest of this section, we discuss equivariancc properties of the ODF expectation. To begin with, 
recall that the translation of a set A C K'' by a G R"* is the set a + A = {a + x:xG A}. Similarly, the 
homothccy of A by a scalar a is given by aA — {ax: x G A}. The homothccy is a dilation for a > 1, 
a contraction for a G [0, 1), and a refiection for a ~ —1. For negative a, the homothccy is a dilation or 
a contraction of the reflection for \a\ < —1 and \a\ G (0,1), respectively. Rigid motion transformations 
g: M'' I— > R'' are given by g(a;) = Fx + a, where F G 0{n) is an orthogonal matrix and a G M''. In W^, these 
transformations are isometrics, i.e. \x — y\ ~ \g{x) — g(j/)|, and form the Euclidean group E{d) with respect 
to the composition. 

Theorem 3.2 (Equivariancc properties). The expected set and the expected boundary are equivariant with 
respect to: 

1. Homothccy, i.e. for a fixed scalar Oi ^ Q, 

E[aA]=aE[A] and E[aaA] = a E[aA]. (12) 

2. The group of rigid motions, i.e. for any g G E{d), 

E[3A] = 5E[A] and E[g 5A] = gE[5A]. (13) 

Proof. ^ We begin by establishing the relation between the ODFs of a set and its homothccy. In particular, 
daA{x) = iTiiy^aA\y ~ x\ = infj,g^ |ay — xj = \a\dA{x/a). Similarly, di^aAyix) ~ \oL\dA':{x/a), so that 
baA{x) = \a\bA{x/a). 

For any y G E[aA], E[6qa(j/)] < 0. Writing y ^ ax for some x, we obtain > E[6qa(j/)] = E[baA{c(x)] — 
\a\ E[6A(a;)]- Thus, E[&A(a;)] < and y ax with x G E[A]. Hence, y G aE[A] and E[q:A] C aE[A]. 

For any y G aE[A], y ~ ax for some x G E[A]. Then E[baAiy)] = ^[\a\bA{y /a)] ~ |a| E[6A(a;)] < 0. 
Thus, y G E[aA] and E[aA] D aE[A]. 

The equivariancc property for the boundary is proved similarly and we omit the details. 

([2]) A transformation g(.T) = r.T + a, for some orthogonal matrix F G 0{d) and a G M.'^, is a bijection 
and g^^{x) ^ T^^{x - 6) = F'^(a; - b). For sets gA = {g(a;) : x G A} and (gA)'= = {g(a;) : x G A""}, we 
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obtain dgA{x) = inf^gg^ \y — x\ = iidy^A jTy + b — x\ = infyeA \y — T'^{x — b)\ = ^^(g ^{x))- Similarly, 
disAyix) = dA<=(g~^(a-)) and bgA = bA{g~^{x)). 

Recall that gE[A] = {g(a;): x e E[A]} and E[gA] = {x: E[6gA(a;)] < 0}. Hence, for any y G gE[A], 
there exists x e E[A] such that y ^ g{x). From here, E[fogA(y)] = E[bA{g^^iTx + b))] = E[bA{x)] < 0. 
Hence, y S E[gA] and E[gA] D gE[A]. 

On the other hand, any y G E[gA] can be written as y = r.T + b, where x = g^^(y). It then follows that 
> E[&gA(y)] = E[6A(g"^(g(a;)))] = E[6A(a;)]. Thus, any y G E[gA] has a form y = Tx + b, where x G E[A]. 
Hence, E[gA] C gE[A] and E[gA] = gE[A]. □ 

Corollary 3.2. Note that any translation in can be described as an isometry with an identity matrix T . 
Thus, it immediately follows that the expected set and the expected boundary are translation- equivariant, i.e. 
for some fixed a G M'', 



E[a + A] = a + E[A] and E[9(a + A)] = a + E[aA] 



(14) 



Theorem 3.3 (Convexity preservation). Suppose that an r.c.s. A is convex a.s.. Then the expectation 
of A is also convex. 

Proof. The a.s. convexity of A implies that cl A is convex and ^a = ^ci a is a convex function a.s. 



(Delfour and Zolesio 



2001 



Theorem 7.1). Hence, for all a G [0,1], bA{ax+{l — a)y) < abAix) + {l — a)bAiy) 
a.s.. Taking expectations, we obtain that E[6a] is convex, and so is E[A]. □ 



We refer to Theorems 13.21 and 13.31 as the shape-preservation properties of the ODF expectation. These 
qualities are particularly desirable in shape and image analysis. In contrast, given a nonatomic probability 
space, the selection expectation of an r.c.s. is convex and coincides with the expec tation of its conve x hull. 
Notably, the convexification of Eg [A] holds even for nonconvex deterministic sets (jMolchanovl . 120051 ) . The 
distance- average ex pectation preserves the convexity , but only if the set, its representative function and the 
window are convex ( Baddelev and Molchanovl . 119981 ) . 



Example 3.2 (Set and its boundary). Consider an r.c.s. A such that A = {0, 1} or A = [0, 1] in IR with 
probabilty p and I — p, respectively. Note that the boundary dA = {0, 1} a.s. in R, so essentially, we are 
observing either the set or its boundary. The expected ODF is given by 



E[6A(a 



= < 



-x, 

(2p- 1)2;, 

(2p-l)(l 
X — 1, 



x), 



X <0 

< a; < 0.5 
0.5 < a; < 1 
X > 1. 



For p = 0.5, E[A] = E[5A] = [0, 1]. For p < 0.5, E[A] = [0, 1] and E[5A] = {0, 1}. For p > 0.5, E[A] = 
E[aA] = {0, 1}. Note that for p ^ 0.5, 9A {0, 1} a.s. However, E[aA] ^ {0, 1}, but rather ^A C E[9A], 
as in Proposition 13. H [5|) . 

For p = 0.5, the distance-average expectation is Eda[ A] = [—1/12, 1/6] U [5/6, 13/12], for m uniform, and 
Eda[A] = {0, 1}, for m = with sufficiently large W (jBaddelev and Molchanovl 1 19981) . 
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More generally, consider an r.c.s. A such that, for some deterministic B E A = B with probability p 
and A = dB, otherwise. From 

bsix), X e B", 
bdsix) = <( 0, xedB, 
-bsix), X e int B, 



it follows that 



Hence, 



bsix), X e B", 

E[6a(x)] = <( 0, xedB, 
(1 - 2p)6b(.t), a; e int B. 



for p ~ 


0.5, 


E[A] 


= E[dA] = B, 


for p < 


0.5, 


E[A] 


= B, E[dA] = dB, 


for p > 


0.5, 


E[A] 


= E[c)A] = dB. 



Example 3.3 (Ball with random radius). Consider a random closed ball A in M'* with a fixed center xq 
and a random radius 8. Recall that the ODF of A is given by 5a (x) = \x — xq 1-9 for all X e M^. Hence, 
E[6A(a:^;0)] = \x — xo\ — E[Q] and the expected set is a closed ball centered at xq with radius E[8]. The 
expected boundary E[9A] is a sphere centered at Xq with radius E[8]. 

Example 3.4 (Flashing discs). Let be a Bernoulli random variable with parameter p. Consider an r.c.s. 
A C given by a disc of radius r centered at the origin, with probability p, and centered at a point a G , 
otherwise. The ODF of A is given by bA{x) = |x|I(0 = 1) + |a; — a|I(8 = 0) — r, where I(-) is an indicator 
function. Taking the expectation with respect to 8, we obtain E[6A(a;)] = p\x\ + (1 — p)\x — a| — r. The 
isocontours of the expected ODF are Cartesian ovals with foci {0, a}. 

Figure [2] shows the expected ODF with the expected boundary superimposed in white, for p ~ 0.8 and 
three different values of a. The observed discs A are shown in gray. When the realizations of A do not 
intersect, the expected boundary is nearly a circle contained in the leftmost disc (Figure [51 left). As the 
distance between the centers decreases, the expected set expands. In general, as p ^ 1, the expected ODF 
approaches that of a disc centered at the origin with radius r. Alternatively, as p 1/2, the isocontours 
approach an elliptical shape. 

For p = 0.5, the c -level isocontour of the expected ODF is an ellipse with foci (0, 0) and a and semimajor 
axis r + c. From E[6A(a::)] > |a|/2 — r, we obtain that if \a\ > 2r, then E[A] = 0; if \a\ ~ 2r, then 
E[A] = {a/2}; if \a\ > 2r, then E[A] is the closure of an ellipse with foci {(0, 0), a} and a semimajor axis r. 
Figure [3] shows the expected boundaries (white) overlayed on the expected ODF, for discs with r = 1 and 
p = 0.5. 

In comparison, the selection expectation Es[A] is always a disc centered at a/2 with radius r, irrespective 
of p. The Vorob'ev expectation Ev[A] for p ^ 0.5 is given by the more likely of the two discs. For p = 0.5, 
Ev[A] is the union of the two discs, provided the cardinality of the intersection does not exceed one, otherwise 
Ev[A] is the intersection of the two discs. 
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Figure 2: The expected boundaries (white) for A with p = 0.8, r = 1 and \a\ = 3,2, 1.5 (from left to right), 
superimposed on the gray-scale image of the expected ODF. Realizations of 3A are shown in white. 




Figure 3: The expected boundaries (white) for random discs with p = 0.5, r = 1 and \a\ = 3,2,1.5 (from 
left to right), superimposed on the gray-scale image of the expected ODF. Realized boundaries are shown in 
gray. 
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Example 3.5 (Random half-plane). Let be a random variable and consider a random half- plane in 
given by A = {x: xi < 8} with dA = {x : xi ~ Q}. The ODF of A is bA{x) = xi — Q and the expectation 
has a form E[&A(a;;G)] = xi — E[8]. From here, the expected set is a half-plane intersecting the horizontal 
axis at point E[e] and E[dA] = {x : xi = E[e]}. In this example, Es[A] = Ev[A] = E[A]. 

Example 3.6 (Random upper half-plane). Consider a random upper half-plane in given by A = 
{x: X2 > XitanG} with dA = {x: X2 — XitanO}. The ODF of A is bAix) — |x| sin (0 — tj), where u; — 
arcsin(a::2|x|~i). For 9 ~ \J[a,b] the expected ODF is E[6A(e) (2;)] = 2{b-a)~^ sm{{b - a)/2)|.T| sin ((a -I- b)/2 - w). 
From here, the expected set is a homothecy with coefficient a = {b ~ a)~-'^(2sin (6 — a)/2))^^ of an upper 
half-plane with boundary angle (a + b)/2. Hence, E[A] is an upper half-plane with boundary given by the 
line that passes through the origin and makes angle {a + b)/2 with the positive horizontal axis. Again, we 
have Es[A] = Ev[A] = E[A]. 

4 Separable Oriented Distance Functions 

Consider a parametric closed set, whose geometry depends on the parameter 9 qMP. We call 9 the generating 
parameter and write A(9) for a closed set generated by 9. For example, a disc with a center 7 and a radius 
p has generating parameters (7,p)- We now extend the notion of parametric sets to random closed sets. 
Specifically, let 8 be a p-dimensional {p > 1) random variable defined on the probability space (il,^, P) 
with E |0| < 00. Consider a random closed set given by a mapping from to J^, so that the geometry of 
the set depends on the realized value of 0. Similarly, we call the generating parameter of a r.c.s. A(0). 
For example, a disc with a random center C and a random radius R has generating parameters = (C, R), 
whilst a disc with a fixed center and a random radius R has a generating parameter = i?. The ODF b^^x) 
of A(0) is a function of random variable and so, is a random quantity. We again assume that bA.{x) is 
integrable for some xq G T> and hence is well-defined. 

Definition 4.1 (Separable ODF). The ODF of set A ~ A{9) with parameter 9 is separable, if it has a form 
bA{x;9) = h'^{x)g{9), where h{x) ~ {hi{x), . . . ,hK{x))'^ and g{9) = (gi(0), . . . ,gK{9))'^ for some functions 
hi{x),gi{9), i = l,...,K. 

Definition 4.2 (Separable r.c.s). An r.c.s. A is separable, if there exists a random variable such that 
6a (a;; 0) is separable a.s.. That is 6a (a;; 0) = h'^{x)g{Q) a.s., with h{x) = [hi{x), . . . jhxix))^ and 
g{9) = {gi{9), . . . , gK{9))^ for some functions hi{x),gi{9), i = 1, . . . ,K. 

The following theorem establishes the connection between the expectation of a set and its generating 
parameter. 

Theorem 4.1. Let & be a generating parameter of a separable r.c.s. A(0) and assume that E[|gi(0)|] < 00 
for every i = 1, . . . , K . Then, 

1. If g{9) is convex in every argument, then E[A(0)] C A(E[0]). 

2. If g{9) is an affine function of 9 in every component gi{9) a.s., thenE[A{Q)] — A(E[Q]). 
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Proof. (P) It follows from Jensen's inequality that E[bA{x;e)] ^ [x)Y\g{&)\ > h{x)'^ g(E[Q]). The right- 
hand side of the equation is the ODF of a set with parameter E[8]. For any x G E[A(8)], we have 
h'^{x)g{E[e]) < E[6a(x; 9)] < and x £ cl A{E[e]) = A{E[e]). Thus, E[A(e)] c A(E[e]). 

^ From Jensen's inequality, E[6A(a;;0)] = h"^ {x)g(E[Q]) iff the function gi{0) is afhne in 9 for every 
i = 1, . . . ,K. The expression on the right-hand site is the ODF of a set with parameter E[0], hence the 
statement ([2]) of the theorem follows. 

Note that although the linearity of g{6) is a necessary and sufficient condition for the form of E[bA.{x; 0)], 
it is only a sufficient condition for the separability of E[A(0)]. This is because the expected set depends on 
E[bA{x; Q)] up to a constant, so that for fixed c > 0,E[bA.{x; 8)] and cE[bA.{x; 8)] define the same expected 
set. □ 

Example 4.1. In Example 1 3. 3 i the ODF of a closed ball with random radius bA{x) = \x\ — 8 is separable 
with h{x) = {\x\, —1)"^ and g{9) = (1, 6)'^ , and the expected set is a closed ball with radius E[6]. 

In Example 13.41 the ODF of flashing discs is not separable, and the geometry of the expected set differs 
from that of the observed sets. However, note that both the random set and its expectation can be described 
as Cartesian ovals. 

In Example 13.51 the ODF of a half-plane is separable with h{x) = (xi, —1)"^ and g{9) = (1, 9)^ . Hence, 
the expected set is a half-plane with parameter E[6]. 

In Example 13.61 the ODF of an upper half-plane can be written as 6a (x) = xi sin 8 — 2:2 cos 8, and so is 
separable with h{x) = {xi, —X2Y' g{9) = {sin 9 , cos 9)^^ . The expected set is an upper half-plane whose 
boundary makes angle E[8] with the positive horizontal axis. This example illustrates the second statement 
of Theorem 14. 11 in that E[A(e)] = A(EQ), even for nonlinear g. 

Remark 4.1. Note that if g is invertible in every component, we can rcparametrize the expected ODF 
as E[&A(ic;8)] = h^{x)g{ri), where rji = g~^{Egi{Q)),i = \,...,K. Hence, the expected set and the 
expected boundary have the same geometric form in terms of parameters r\i. The relation between r\i 
and the moments of 6 is established on a case-by-case basis. In particular, for the ODFs of the form 
^k{x) ~ h{x) + g{Q) for some function h{x) : R*^ ^ M and an invertible function g{9) : W M., we obtain 
that E[6A(a;;8)] = h{x)+E[g{e)] = Hx) + g{g-^E[g{e)])) and E[A(8)] - A{g~^ E[g{e)]). 

Remark 4.2. In digital processing, data are stored as an array of pixels and so the ODF of an object in 
the image is discretised. This discretization implies that for realizations of an r.c.s. A at times t ~ 1,2, . . ., 
the ODFs of At can be viewed as a random sample of arrays in the time domain. For a separable ODF, the 
space-time covariance matrix has a form 

cov(6a, (x), 6a, (y)) = h^{x) cov(g(et), g{Qr))h{y), 

for lattice vertices x, y and time points t, r. This expression shows that the covariance of a separable 
ODF is a product of a deterministic spatial component and the covariance of the stochastic time compo- 
nent. In particular, for an ODF of the form bA.{x) = h{x) + g{Q), the space-time covariance is given by 
cov{bAt{x),bA^{y)) = cov{g{Qt),g{'c>T)), and hence the dependence between the sets can be inferred from 
their generating parameters. 
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5 Further Developments and Applications 



5.1 The Sample Mean Set 

Let Ai, . . . , Am be the observed realizations of an r.c.s. A. For every i, we assume that A^ C cl I? and 
dAi ^ and denote by bi{x) the ODF of A^ at point x. Let bm{x) = Y^iiLi be the sample mean 
ODF at point x. 



Proposition 5.1. The sample mean ODF bm{x) is uniformly Lipshitz with constant one. 
Proof. The proof follows immediately from the Lipshitz continuity of biix). 



□ 



Definition 5.1 (Sample mean set). The sample mean set A^ and the sample mean boundary dA-m are 
given by the zero-level set and by the zero-level isocontour of the empirical ODF, respectively, i.e. 



{x : bm{x) < 0}, dAjn = {x : b,n{x) = 0}. 



(15) 



The sample mean set defined in ([TS|) is a random closed set, since {A„i n if 7^ 0} = {mix^K bm{x) < 0} is 
measurable. 



Sample size, m 



Sample size, m 



Figure 4: (Left) The ratio of the radii of A^ to that of the expected set E[A] in Example 15.11 for 8 ~ 
U[0.8, 1.2] (black) and 8 ~ U[0.5, 1.5] (gray). (Right) The difference between the boundary angles of A™ 
and E[A] in ExamplcOwith 6 ~ U[7r/8, Stt/S] (black) and 8 - U[0,7r/2] (gray). The median (solid) and 
the 25-tli and 75-th percentile (dashed) values are based on 1000 simulations for each m. 



Example 5.1 (Disc with random radius). Recall that for a closed ball A with a random radius 8, E[A(8)] = 
A(E[6]) and E[3A] = 9A(E[6]). For a sample of closed balls Ai, . . . , A™ with random radii 81, . . . , 8m, 
the emprirical ODF is bm{x) = |a;| — 8™ and hence A™ — A(8m) and dAm — 9A(8m). 

For an independent random sample of radii values 8i,...,8,„, we construct m random discs in 
centered at the origin. To compare the empirical boundary estimator with the theoretical average, we 
compute the ratio of the radii of E[A] to A(8m); see FigurelU left. The radii values were sampled from the 
uniform distribution U[0.8, 1.2] (black) and U[0.5, 1.5] (gray). The median values of the ratios (solid curves) 
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and the 25-th and 75-th percentiles (dashed curves) are based on 1000 samples for each sample size m. The 
figure shows that the accuracy of the empirical estimate is higher for bigger sample sizes and for smaller 
variance values of the generating parameter. 



Example 5.2 (Random upper half-plane). Consider an r.c.s. in Example 13.61 given by A(0) = {a- 6 
M.'^ : X2 > a;itan0}, where Q is the angle between the boundary of the plane and the positive horizontal 
axis. Recall that for 6 - U[a,6], E[A(e)] = A(E[e]) and E[aA] = dA{E[e]). 

For an independent random sample of the angle values 0i, . . . , 8„i, we construct m upper half-planes. 
The empirical mean boundary is estimated by the zero- level isocontour of the empirical mean ODE. Figure |4] 
(right) shows the difference between the boundary angles for the empirical mean set and that of the expected 
set for 9 ^ U[7r/8, Stt/S] (black) and Q ^ U[0, 7r/2] (gray). The median values of the ratios (solid curves) 
and the 25-th and 75-th percentiles (dashed curves) are based on 1000 replicates for each value of m. As in 
the previous case, the boundary estimate is more accurate for bigger sample sizes and for smaller variance 
values of O. 

The empirical mean set and the empirical boundary provide a recipe for practical construction of the set 
and boundary estimators. Examples 15.11 and 15.21 suggest the consistency of the boundary estimator, based 
on the empirical ODE. Note that both the expected set E[A] and the empirical mean set A,„ are described 
as level sets of the continuous function 6a and hence , unde r certain conditions, the coi isistency of Am. a s 
an estimator of E[A] can be inferred from (jMolchanovl . Il998l . Theorem 2). The result of (jMolchanovl Il998f ). 
however, does not apply to the boundary es timator dAm.. In our upcoming paper, we study the consistency 
of the boundary estimator in greater detail I Jankowski and Stanberrvl (|in preparatiorJ l . 



5.2 Loss Functions 

Assume that the observed data A G follows a distribution Fa with mean parameter 2t = E[A], 21 e JT. We 
call 2t the parameter set. Here, data A represent random closed sets with nonempty boundary. In practice, 
it is more common to acquire images rather than observe sets per se, so the observed A can be seen as a 
result of the first-step analysis. The goal is to estimate the parameter set 21 given the observed data A. 

If set 2t is parametric with parameter 0, we can construct an estimator of 21 using conventional loss 
functions for 9. Eor example, if A is a collection of deformed balls, the distribution F can be parametrized 
by the center and radius of a ball. Depending on the problem, we might be interested in testing the hypothesis 
about either the location or the size of the ball, or both. Then the optimal estimator 6 is determined using 
standard loss functions, e.g. the indicator loss £{6, a) = I{a = 6), the absolute error loss i{d,a) = \a — 9\ or 
the squared error loss £{0,a) = (a — 0)^. However, often, the sets of interest cannot be parametrized and 



at best can be described as "blobs" , in which case conventional loss functions are not applicable. Here, 


relay loss functions for Bayesian inference about sets introduced in 


Stanberrv and Besae 


in oreoaration) 



the frequentist framework. 



To begin with, we consider a loss function based on representation of a set A given by its characteristic 
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function Xa{x) = 1a{x). The space of characteristic functions is complete with metric 
d{A,B) = ||X.4-X.,I|l. - 



1/9 

\XA^)-XsixW\idx)] . (16) 



Here and later, the integration domain V is either the entire image or its subset. Let X be an estimator of 
Xa- Given the metric d, we define an Li loss function as 



h{X^,X)^\\X^-X\\L^= f \X~Ax)-X\X{dx). 



(17) 



Because the estimator X is the characteristic function that minimizes the expected loss, the estimator 2t of 
21 is immediately available. Strictly speaking, d{A, B) distinguishes between sets up to a set of Lcbesgue 
measure zero, so we assume that the parameter set has no punctures and hence 21 is uniquely determined. 
Note that the equation (fT7|) is the measure of the symmetric difference, since 

4(X=i,X) = A(2lAi) = A((2l\2i)U(2i\2l)) (18) 

and, hence, ti{Xn,X) is akin to the absolute error loss in conventional settings. The set tha t minimizes 
E[A(2 tA2i)] is the Vorob'ev median, i.e. is the median level set of the coverage function ([3]); (see 



Molchanov 



2005 



p. 178). The Vorob'ev median is a set-analog of the ordinary median in the sense that, as the latter, 
the former minimizes the absolute error loss. 

Moreover, the Vorob'ev median is optimal with respect to a much broader class of Lq loss functions, 

^g(X«,X) = llXa -XIIl,, 1<<z<^. (19) 
Thi s is because Xn takes values in {0, 1} and the induced Lq topologies (fT9|) are all equivalent for 1 < g < oo; 



(see lPelfour and ZolesioL I2OO1I Theorem 2.2). 

Instead of optimizing the expected Lq loss over the set of characteristic functions, consider minimizing 



over a larger set of functions with values in [0, 1]. For the squared L2 loss 

h{X^,X) = \\X^-X\\l,, (20) 

the pointwise estimator X{x) G [0, 1] is the coverage function ([3]). The estimator X{x) is a global minimum in 
L2 since E[||Xg — E[X2j|||^] < E[||X2i — XH^^]. Note that X is essentially a gray-scale image with intensities 
reflecting the probability of a pixel being in 21. 

The coverage function, however, is not necessarily a characteristic function, unless A is deterministic 
a.s.. Consider a set- valued estimator given by an excursion set of X{x) so that determining the optimal 
estimator 21 now reduces to choosing an appropriate threshold u for the excursion set ([3]). The Vorob'ev 
criterion ([5]) gives an estimator 21, which is optimal with respect to Lebe sgue measu r e in t he sense that 



E[A(2lA2l)] < E[A(2tA2t)] for all measurable sets with measure E[A(A)] (see lMolchanovll2005l . Theorem 2.3, 
pg. 177). From p8|l and the equivalence of the ig-topologies on the space of characteristic functions, it 
follows that the Vorob'ev expectation minimizes the expected loss for any g G [1, 00) over all of the sets 
with Lebesgue measure E[A(A)]. 
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In general, given the representative functions /a and of sets A and B, the discrepancy between A and 
B can be taken in terms of some (pseudo-) metric 



m{A,B)^mw{fAjB), (21) 

defined on a compact window W C P. Further, we suppress the subscript W whenever possible. Given the 
metric m{A, B) = \ \fA ~ /b||l2, the global estimator based on the squared L2 loss 

^2(/a,/) = ll/2i-/llL (22) 

is the average function E[/a]- Hence, the distance-average expectation ([6]) is an optimal estimator of 21 with 
respect to (P^. among all of the level sets of E[/a]. 

Note that although the distance-average expectation arises naturally for the L2 loss, it was originally 



defined for a generic (pseudo-) metric m(-, •) (jBaddelev and Molchanovl . Il998f) . Thus, in principle, we can 



construct an estimator as the level set of the expected representation for arbitrary choice of m. For that, 
we must specify the representative function /, select an appropriate (pseudo-)metric m, and choose the loss 
function £(/a, /) based on m. The threshold u in ([6]) then determines the optimal estimator with respect to 
£{f<2i, /) among all of the level sets of E[/A(a;)]. 

We now consider a representation of a set given by its ODF. Recall that the space of ODFs of sets with 
non-empty boundary is complete with metric d{A, B) = — ^aHl,- For the loss function 

^2(6a,^) = = / \b<^{x)-b{x)\^X{dx), (23) 



let b be the estimator of b<^ that minimizes the expectation of ((23|) . Because the ODF gives a unique (up to 
the boundary) representation of a set, the estimator 21 is given by the zero-level set of 6. However, optimizing 
E[^2(^a, b)] over the set of ODFs is nontrivial, so instead we consider unrestricted minimization over the set 
of functions defined on T>. In this case, the pointwisc estimator is the expected ODF, which is also a global 
minimum. Because the estimator E[6A(a;)] is not necessarily an ODF itself, it does not uniquely determine 
the set. Hence, we take the set-valued estimator 2t of 21 to be the zero-level set of E[&A(a;)], so that 2t is the 
expected set as defined in p.ip . The estimator of the boundary 921 is given by the zero- level isocontour of 
E[bAix)]. 

5.3 Image Averaging 



In thi s section, we consider the example of image averaging originally discussed in lBaddelev and Molchanov 



(jl998l ). In image averaging, which relates to Bayesian image classification and reconstruction, the goal is is 
to determine an average object or a typical shape from the collection of images of the same scene or objects 
of the same type. 



To begin with, we quickly outline the image sampling procedure and refer for more details to lBaddelev and Molchanov 



(|l998l ). Figure [5] (left) shows the true binary image, /, of the scanned and thresholded newspaper fragment. 
After adding Gaussian noise, 15 independent realizations from the posterior were obtained using a Gibbs 
sampler with true noise parameters in the likelihood and an Ising prior. In the current notations, I? = W = /, 
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21 is the true newspaper text a nd Ai , . . . , A15 a re the 15 independent reconstruction of 21. The data for this 
example was downloaded from lBaddelev ( 2009 ). 




Figure 5: The true binary image of the newspaper fragment. 

Figure [5] shows the ODF estimator ^odf (left) and the distance- average estimator 2 l^)/^ (right) with 
/a(-) = bAi') and m(-,-) = L2{-,-); see also Figure 6 in iBaddelev and Molchanovl (|l998l ). The estimator 
2I0DF appears less noisy and more accurate as compared to ^da, where the letters look overinflated. 
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Figure 6: The ODF estimator 21odf (left) and the distance-average estimator 21d^ (right). 



To compare the results of the reconstruction, it is instructive to look at the residual images for ^odf 
and 'QioA in Figure[7l left and right, respectively. Note that black (resp. gray) pixels in the figure correspond 
to the set 21 \ 2t (resp. 2t \ 21). The residual image for ^da shows a clear spatial pattern, indicating that 
'Qi-DA tends to overestimate the set. In turn, the residuals of 21odf show little spatial clustering. In addition, 
2I0DF C ^DAi that is the distance-average estimator is inclusive of the ODF one, with discrepancies observed 
along the boundary and in the background. 

To compare the quality of the estimators, it is common to use the fraction of misclassified pixels. Here, the 
error is 4.12% for 'QioDF and 10.68% for pipA- However, a small misclassification error does not necessarily 
imply similar looking images (jBaddelevl . Il992f ). therefore, we compare the two estimators by computing an 
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Figure 7: The residual images of the discrepancy between the estimated set and the truth for ^odf (left) 
and ^DA (right). Black (resp. gray) pixels mark set 2t \ 21 (resp. 2l\2t). 



Lo distance to the truth 



d(2l,2t) 



(/ Hix) 

JV 



6a(x)p dxf'\ 



which is compatible with the distance- average expectation. Calculating the distances for the two estimators, 
we obtain (i(2l, 21ol>f) = 1-78 and rf(2l, SIda) = 2.66, which shows that the ODF estimator outperforms the 
distance-average reconstruction. 

Note that the reconstruction results using ODFs do not change if we reverse the color of pixels, as 
opposed to the dist ance-average estimator wi t h //; ( ■) = rfyi(-)j where the swapping affects the quality of the 
reconstruction (see 



Baddelev and Molchanov 



19981 Figures 4 and 5). 



5.4 Boundary Reconstruction in Noisy Images 

In image analysis, it is often required to reconstruct the boundary of an object in a noisy image. Low level 
methods such as edge detection and intensity thresholding strongly d epend o n image characteristics and 



1986 



Gonzalez and Woods 



requi re a careful tuning of parameters to produce satisfactory results ([Canny 
20021 ). A number of methods for boundary reconstruction have been proposed in the literature. In active con- 
tour models, the boundary is represented by a parametric or a free- form curve and the recons truction is driven 



by local image force s , with smoothness control l ed via user- imposed regularity constrain ts ( Blake and Isard 



1998 




Brigger et al. 




2000 


■ 


El-Baz et al. 




2005 




Flickner et al. 




1994 


Kass et al. 




1988 


Rueckert and Bureei 




1997 




Tauber et al 




2004 




TerzoDOulos et al.. 


1987 




Yuille et al. 



Menet et al 



199C; 



19921 ). The use of 



local properties of the image implies strong dependence of the boundary estimate on initialization of the 
modeling curve, a contrast-to-noise ratio of the image and convexity of the object. Various extensions were 



proposed in the literature to improve the performance of acti ve contours (jCoughlan et al 



2005 



Wang et al 



1996 : 



Cohei 



1991 



Xu and Prince 



19981) 



200C; 



El-Baz et al 



Many Bayesian m ethods for boundary reconstruction use templates which rep r esent local or global ge- 



ometry of an object ( Grenander 



1976 



Tielineland and Besag 



1998; 



Amit et al 



1991 



Grenander 



1996' 



Jain et al 



1996 : 



Grenander and Keenan . 



1993 : 



Grenander 



1996 



Jain et al 



1998 : 



Chalmond 



2003 : 



Geman et al 
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19901 ). The stochastic method in lQian et alj (|l996f ) models the boundary by closed polygons. The method 
was then extended to allow for the dynamic es timation of the number of polygon vertices, thus balancing 
the complexity and the flexibility of the model (jPievatolo and Greenl . 119981 ). However, the reconstruction is 
given by a gray-scale image, where intensities reflect the posterior mean state of pixels, so that constructing 
a bou ndary estimator requi r es additional po st-processing. 



In 



Stanberrv and Besagl (|in preparatiorJ l . the authors proposed a new method to reconstruct a smooth 



connected boundary of an object in a noisy image using B-spline curves. The method is Bayesian and uses a 
Markov chain Monte Carlo algorithm to obtain curve samples from the posterior. Using the squared _L2-loss 
based on ODFs, the posterior estimator is given by the posterior expected boundary as defined in (|10p . 
The simulation study showed that the ODF estimator is more accurate as compared to the distance-average 
reconstruction and Vorob'ev sets. 



5.5 Implementation 

In addition to its appealing theoretic al propertie s , the PDF expectation can be efficiently computed, usin g 



algorithms for the distance function (jBreu et al 



1995 



Freidman et al 



1977 



Rosenfeld and Pfaltz 



19661). 



The ODF is computed by combining the distance function for the original image and that of the inverted 
image, i.e. the black and white colors are interchanged. The examples and results reported in this paper 
were coded in MATLAB (R2008A, The MathWorks), where the distance function to a set is computed using 
the bwdist command. 

In comparison to the distance-average and Vorob'ev expectations, the ODF average is more efficient, 
because it requires no optimization since the threshold level is fixed at zero. This difference is particularly 
prominent in Bayesian reconstruction, where the estimator is based on thousands of samples from the 
posterior. In addition, the storage requirements for the ODF reconstruction in Bayesian framework are 
miniscule, as we only need to update the average after every sweep. 



6 Discussion 

In this paper, we present new definitions of the expected set and the expected boundary, using the ODF 
representation of a set. In conventional settings, where the parameter of interest is a scalar or a vector, it is 
straightforward to construct an estimator of the mean from the observed data. However, statistical inference 
about sets is nontrivial. When dealing with sets, one has to determine the features that are important to 
emphasize and develop inference methods based on them. These features may include the location, the size, 
or the orientation of a set. If a set has an analytic representation, the most natural approach to set inference 
is via inference about its parameters. This, however, is not applicable for sets with arbitrary geometry. 
A number of existing definitions of the expected set are based on the linearization idea outlined in 



Molchanovl (|2005l ). where inference is made on the space of functions representing a set. The resulting 
estimator necessarily depends on the type of linearization and the choice of the representative function. In 
this paper, we use an approach similar in spirit to the linearisation idea and define the expectation based 
on the ODF representation of a set. However, we forgo the optimality criterion and simply choose zero as 
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a threshold. Whilst the choice might appear simplistic, it gives an expectation with attractive theoretical 
traits, including cqui variance, convexity-preservation and inclusion properties. 

The definition is particularly appealing for problems in image analysis. In particular, its denoising prop- 
erty implies that random specks will be averaged out unless observed with probability one. In comparison, 
the distance-average expectation is not empty by construction, and therefore identifies false positives in any 
noisy image. Furthermore, the distance-average reconstruction strongly depends on the choice of pseudo- 
metric, representative function, restriction window and any parameters of the above. In turn, the selection 
expectation depends on the structure of the probability space, whilst the Lebesgue measure criterion of the 
Vorob'ev expectation is not well suited for image analysis, although, the optimality of the Vorob'ev median 
is rather attractive. 

In our upcoming paper, we study the asymptoti c properties of the PDF estimators for the boundary and 



describe a method for constructing confidence sets (jjankowski and Stanberrvl . lin preparation! ). 
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